Import reuireed modules


In [7]:
import sys
sys.path.append('../')

import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline

from porousmedialab.column import Column

Initiate the model and domain


In [8]:
tend =3
dt = 0.0025
dx = 0.2
L = 10
phi=0.99
ftc = Column(L, dx, tend, dt, w=5)

Add species


In [10]:
C_init=0.0000001/phi/3
ftc.add_species(theta=phi, name='H2CO3', D=50, init_conc=C_init, bc_top_value=0, bc_top_type='flux', bc_bot_value=0, bc_bot_type='flux')
ftc.add_species(theta=phi, name='HCO3', D=50, init_conc=C_init, bc_top_value=1e-3, bc_top_type='dirichlet', bc_bot_value=0, bc_bot_type='flux')
ftc.add_species(theta=phi, name='CO3', D=50, init_conc=C_init, bc_top_value=0, bc_top_type='flux', bc_bot_value=0, bc_bot_type='flux')

ftc.add_acid(species=['H2CO3', 'HCO3', 'CO3'], pKa=[3.6, 10.32])


ftc.add_species(theta=phi, name='H2SO4', D=50, init_conc=C_init,  bc_top_value=0, bc_top_type='flux', bc_bot_value=0, bc_bot_type='flux')
ftc.add_species(theta=phi, name='HSO4',  D=50, init_conc=C_init,  bc_top_value=0, bc_top_type='flux', bc_bot_value=0, bc_bot_type='flux')
ftc.add_species(theta=phi, name='SO4',  D=50, init_conc=C_init,  bc_top_value=1e-6, bc_top_type='dirichlet', bc_bot_value=0, bc_bot_type='flux')

ftc.add_acid(species=['H2SO4', 'HSO4', 'SO4'], pKa=[-10, 1.99])

ftc.add_species(theta=phi, name='Ca', D=50, init_conc=1e-9, bc_top_value=1e-2, bc_top_type='dirichlet', bc_bot_value=0, bc_bot_type='flux')
ftc.add_species(theta=1-phi, name='CaCO3', D=5, init_conc=0., bc_top_value=0.0, bc_top_type='flux', bc_bot_value=0.00, bc_bot_type='flux', w=0.)
ftc.add_species(theta=1-phi, name='CaSO4', D=5, init_conc=0., bc_top_value=0.0, bc_top_type='flux', bc_bot_value=0.00, bc_bot_type='flux', w=0.)


ftc.add_ion(name='Ca', charge=2)

Add rate constants


In [11]:
ftc.constants['Ks_CaCO3'] = 3.3e-9
ftc.constants['Ks_CaSO4'] = 10**-4.58
ftc.constants['k_pre'] = 1e-4
ftc.constants['k_dis'] = 1e-3
ftc.constants['CF'] = (1-phi)/phi  ### conversion factor

Rates


In [12]:
ftc.rates['R_pre_CaCO3'] = 'k_pre * (Ca*CO3/Ks_CaCO3-1)'
ftc.rates['R_pre_CaSO4'] = 'k_pre * (Ca*SO4/Ks_CaSO4-1)'
ftc.rates['R_dis_CaCO3'] = 'k_dis * CaCO3 * (1 - Ca*CO3/Ks_CaCO3)'
ftc.rates['R_dis_CaSO4'] = 'k_dis * CaSO4 * (1 - Ca*SO4/Ks_CaSO4)'

ODE for particular species


In [ ]:
ftc.dcdt['CaCO3'] = 'R_pre_CaCO3/CF - R_dis_CaCO3'
ftc.dcdt['Ca'] = '-R_pre_CaCO3 + R_dis_CaCO3*CF - R_pre_CaSO4 + R_dis_CaSO4*CF'
ftc.dcdt['CO3'] = '-R_pre_CaCO3 + R_dis_CaCO3*CF'
ftc.dcdt['CaSO4'] = 'R_pre_CaSO4/CF - R_dis_CaSO4'
ftc.dcdt['SO4'] = '-R_pre_CaSO4 + R_dis_CaSO4*CF'

In [ ]:
ftc.solve()


Simulation started:
	 2018-09-14 14:13:12


Estimated time of the code execution:
	 0h:00m:54s
Will finish approx.:
	 2018-09-14 14:14:11

In [ ]:
ftc.plot_contourplots()

Saturation indeces


In [ ]:
ftc.plot_saturation_index('Ca','CO3','Ks_CaCO3')

In [ ]:
ftc.plot_saturation_index('Ca','SO4','Ks_CaSO4')

In [ ]: